Effect of Gravity and Confinement on Phase Equilibria: 
A Density Matrix Renormalization Approach 
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The effect of gravity on a two dimensional fluid, or Ising magnet, confined between opposing walls 
is analyzed by density matrix renormalization. Gravity restores two phase coexistence up to the 
bulk critical point, in agreement with mean field calculations. A detailed finite size scaling analysis 
of the critical point shift is performed. Density matrix renormalization results are very accurate 
and the technique is promising and best suited to study equilibrium properties of two dimensional 
classical systems in contact with walls or with free surfaces. 
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The thermodynamical properties of confined systems 
have received a lot of attention in the past years [P^0| . 
Their critical behavior is rather different from the bulk 
criticality and has been the subject of extended inves- 
tigations by means of mean field and scaling analysis 
Monte Carlo simulations J5|-|| and exact calcula- 
tions ||[l0]]. The simplest and most studied case is the 
Ising model in a L x M d ^ 1 lattice with M — > oo, i.e. 
confined between two infinite walls separated by a finite 
distance L. Of considerable interest is the situation in 
which magnetic fields (hi and h 2 ) acting on the spins at 
the walls are introduced. 

For parallel surface fields (hi -hi > 0) and finite L 
two phase coexistence is shifted to finite values of a bulk 
magnetic field ft, O, as illustrated in Fig. |l|(a). This phe- 
nomenon is analogous to the capillary condensation for 
a fluid confined between two parallel surfaces, where the 
gas-liquid transition occurs at a lower pressure than in 
the bulk. Finite size scaling predicts that the capil- 
lary critical point [h c (L),T c (L)} scales as: 



T C (L) - T c ~ L~ VT and h c (L) ~ L~ 



(1) 



where T c is the bulk critical temperature and i/t = 1/f, 
jjh = d — (3 /v are the thermal and magnetic exponents, 
respectively (here d is the dimensionality, v and /3 are the 
correlation length and magnetization exponents, respec- 
tively) . 

At fixed T < T c and finite L the scaling to the first 
order line is of type (!]] : 



h (L) ~ 1/L. 



(2) 



While the previous relation has been verified in Monte 
Carlo simulations in d = 2,3 a direct verification 

of the scaling of the capillary critical point has not been 
attempted yet, due to the high computational effort |t]] 
needed to locate accurately [h c (L),T c (L)\. 

The case of opposing surface fields (for simplicity we 
consider hi = —h%) was analyzed in detail by Parry 
and Evans using a Ginzburg-Landau approach. They 



found that two phase coexistence is restricted to tem- 
peratures below the interface derealization temperature 
Td(L) as shown in Fig. 0(a); the surprising result 0] 
is that Td(L) does not scale to the bulk point T c , for 
L — > oo, but to the wetting temperature T w as: 

-1/0S 



T d (L) - T w ~ L~ 



(3) 



Here (3 S is the exponent describing the divergence of the 
thickness of the wetting layer for a semi infinite system: 
I ~ (T w — T)~@ s . We recall also that T w depends on the 
value of the surface field hi and can be far away from the 
bulk critical point (see Fig. g(b)). For T d (L) < T < T c 
there is a single phase with an interface meandering 
freely between the walls. Numerical results from exten- 
sive Monte Carlo simulations [|[D in d = 2, 3 confirm the 
mean field scenario. 
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FIG. 1. (a) Phase diagram of the d-dimensional Ising 
model for a bulk system in the (h,T) plane (solid line). The 
dashed lines are the phase diagrams of confined systems with 
identical (1) and opposing (2) surface fields, (b) Dependence 
of the wetting temperature on the surface field h\ . 

Rogiers and Indekeu Q extended the model with op- 
posing surface fields, including a bulk field which varies 
linearly along the finite direction of the system and mod- 
els the effect of gravity on a confined fluid. They consid- 
ered the following Hamiltonian: 
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where J > 0, S(jj) = ±1 and 1 < i < L labels the 
lattice sites along the finite direction, while j labels the 
remaining d— 1 directions. The parameter g is the equiv- 
alent of the gravitational constant. Using a mean field 
approach, it was found Q that the ordinary scaling to the 
bulk critical point is restored in an extended parameter 
space where g is included. While the mean field analysis 
is correct at dimensions higher than the upper critical di- 
mension (d = 4), its validity in lower dimensions, where 
thermal fluctuations become important, is questionable. 

The aim of this Letter is to study the phase diagram of 
the model described by the Hamiltonian (jij) for anixoo 
strip, beyond the mean field approximation. We test the 
validity of the conclusions of Ref. jl| at the lower critical 
dimension (d = 2) and clarify the mechanism of the crit- 
ical point shift for a system confined between opposing 
walls. We analyze also the exponents governing this shift 
and test the validity of scaling assumptions that have 
been proposed some time ago but, to our knowl- 

edge, never derived directly from a microscopic model. 

The numerical calculations are based on a density ma- 
trix renormalization group (DMRG) approach, a tech- 
nique developed by White |T2|]lg| ] for the study of ground 
state properties of quantum spin chains. Nishino ex- 
tended the DMRG to d — 2 classical systems as a transfer 
matrix renormalization. Transfer matrices are frequently 
used for numerical investigations of L x oo strips; compu- 
tations are restricted to small strip widths since the di- 
mension of the transfer matrix grows exponentially with 
L p!5| . In the DMRG algorithm one constructs effective 
transfer matrices of small dimensions (i.e. numerically 
tractable), but which describe strips of large widths. The 
spin space is truncated in a very efficient manner and the 
numerical accuracy is very good p2|-^4|, even for renor- 
malized matrices of not too large dimensions Jig . More 
details will be presented elsewhere [ fL7| . 

We start the analysis of the phase diagram of the model 
from its ground state properties. At T = 0, hi < J and 
for a strip of width L the ground state is double degen- 
erate (two phase coexistence) with all spins either up or 
down for: 

jjs(hi-J) < 9 < ^(h 1 + J) (5) 

For values of g outside this interval the ground state is 
non degenerate, with an infinitely long straight interface 
at the center of the strip separating a region with all spins 
up from a region with all spins down. 

To distinguish the two phase coexistence region from 
a single interface-like state it is convenient to calculate 
the correlation function between two neighboring spins 



at the center of the strip c L / 2 = ( s (i j) s (i + i At 
T = 0, cl/2 drops from +1 in the two phase coexistence 
region to —1 in the single phase region. For nonzero tem- 
peratures, a sharp distinction between these two regions 
is not possible since true criticality does not occur for 
finite values of L. However a peak in the temperature 
derivative of c^/ 2 can be seen as a pseudocritical point 
where the system goes through a smooth change from 
two phase coexistence to a one phase state, without a 
true phase transition. As the strip width increases these 
finite peaks shift towards the true thermodynamic sin- 
gularities of the bulk system. The analysis of shifts of 
pseudocritical points has been considered already 
for the present model in absence of gravity. 
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FIG. 2. Phase diagram of the model for different strip 
widths in the (g,T) plane for hi — 0.5 and J = 1; the area 
below the curves is the two phase coexistence region. Inset: 
Scaling of T ma x(£) and Td(L) vs 1/L. Error bars are much 
smaller than symbol sizes. 

Figure || shows the phase diagram in the (g,T) plane 
for J = 1, hi = 0.5. Each curve is the phase bound- 
ary between the two phase coexistence (area below the 
curve) and the one phase region, for a specific value of 
L. Only the phase boundaries for L up to 40 are shown 
although the largest size considered is L = 60. For g > 
two phase coexistence is shifted to higher temperatures 
with respect to g = due to a competing effect of surface 
and bulk fields, while at negative g phase coexistence is 
suppressed. As L increases the whole two phase coexis- 
tence region shrinks and shifts towards the g = axis. 
The intersections of the phase boundaries with this axis, 
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indicated by the horizontal arrows in Fig. 0, define the in- 
terface derealization transition temperatures Td(L); the 
phase boundary maxima [g max (L) ,T max (L)] are indicated 
by the vertical arrows in the figure. 
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FIG. 3. Plot of ln[(; m ax(i)] vs InL for different surface 
fields. The dotted lines correspond to a slope —2.875. Er- 
ror bars are smaller than symbol sizes. 

Wc find for Td(L) a scaling in excellent agreement with 
Eq. (||) [||,|| and for T max (L) a scaling of the type: 



7max(£) — T c ~ L 



II T 



(6) 



with the d — 2 Ising exponents /3 S — 1 and j/t = 1. A plot 
of T max (I/) and Td(L) for /ii = 0.5 as function of 1/L is 
shown in the inset of Fig. 0. From an extrapolation of the 
data in the figure we find the following estimates for the 
wetting and bulk critical temperatures: T w — 1.9584(8) 
and T c = 2.272(4) (the exact values are T w = 1.95845 
p8| and T c = 2.269). The same scaling analysis has been 
extended to several values of h\. all data points are in 
excellent agreement with the relations (|3|) and (^j). 

Van Leeuwen and Sengers ]ll[ analyzed the influence 
of gravity in an infinitely extended system; they argued 
that the product of g with a length in the direction of 
the gravitational field should scale as a bulk constant 
field. On a finite strip, this leads to the following scaling 
prediction Q: 



gmax(L) 



(7) 



Both at T = (see inequality (|5j)) and also for finite T 
Jl7t , not too close to the phase boundary maxima, we find 
a scaling towards the bulk first order line of type 1/L 2 , 
in agreement with Van Leeuwen and Sengers' conjecture 
since the corresponding scaling of a constant bulk field is 
1/L (see (|)). 

Figure pTshows a plot of ln[g max (L)] vs ln(L) for four 
values of hi, from 0.1 to 0.99 and J = 1 pi. The dotted 



lines, drawn as a guide to the eye, correspond to a scaling 
of type (0) with the d = 2 Ising exponent 1 + yn = 2.875. 
The data agree with the scaling relation (^) for the largest 
surface field considered (hi — 0.99) for which a linear fit 
of the points for L > 20 yields an exponent 2.86(2) in 
good agreement with the Ising value. At lower surface 
fields and up to the largest size considered (L = 60) the 
exponent deviates from the expected value and increases 
from 2.4 to 2.6 as the surface field increases from hi = 0.1 
to hi — 0.5. The scaling analysis for <? max (£) is some- 
what more problematic than the scaling of T max (L); this 
was to be expected since along the gravitational field di- 
rection there is a scaling to T c and a scaling to the bulk 
first order line (of type 1/L 2 ). The interplay between 
these two scalings may be the cause of the observed shift 
of the exponent of g niax (L) from the Ising value for low 
surface fields, where the asymptotic behavior ((7J) possibly 
sets in for L ^ 60. In other studies of confined systems it 
was found that, in the scaling analysis of finite size data, 
the value of the system width beyond which one has a 
clear asymptotic behavior, may depend strongly on the 
surface field (see for instance Ref. 0). 

Figure [| shows some magnetization profiles for L = 60, 
J = 1, hi = 0.5, T = 2.0 and for different values of g. 
The profile (a) corresponds to a point of the phase dia- 
gram located in the two phase coexistence region, where 
the magnetization is averaged over two phases. This pro- 
file is similar to those calculated exactly by Maciolek and 
Stecki 0,|2fJ for g = and T < T d {L). Profiles (b) and 
(c) correspond to points of the phase diagram in the one 
phase region. Notice the competing effect between grav- 
ity and surface fields in the vicinity of the walls in the 
profile (b). 

Magnetization profiles of interface-like configurations 
can be calculated using a solid-on-solid (SOS) approxi- 
mation pl[. It is assumed that all the spins at the two 
sides of the interface are fixed and take values ±rao; in- 
terfacial configurations with overhangs are neglected. 

SOS magnetization profiles are shown in Fig. |] (dashed 
lines) and are given by jljj : 
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x(i)/ti 



dt e~ l 



(8) 



where x(i) = L[i-(L+l)/2]/(L- 
width, is given by the relation: 



1) and the interface 
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with (To the surface tension. For <ro and mo we take their 
exactly known values at g = 0; this approximation is 
good for low gravity. 

Results from the interface Hamiltonian agree well with 
those calculated with DMRG even at temperatures not 
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too far from T c (in Fig. 4 the temperature is approx- 
imately 10% below T c ). Gravity has the effect of re- 
ducing interface fluctuations, so that configurations with 
overhangs, which are not taken into account in the SOS 
model, have small weights. We stress that (||) is valid in 
the limit £j_ <C L since the effect of the walls has been 
neglected. 
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FIG. 4. (Circles) Magnetization profiles in the two phase 
coexistence region (a) and in the single phase region ((b) and 
(c)) calculated by DMRG (Dashed lines) Profiles calculated 
with a SOS approximation [Eq. |8j. Inset: Detail of the profile 
(b) far from the interfacial region. 

In conclusion, we found that the competing effects of 
surface fields and gravity restore the "ordinary" finite 
size scaling to the bulk critical point, in agreement with 
mean field predictions which thus survive the strong 
thermal fluctuations at the lower critical dimension. It is 
currently believed that the model confined between op- 
posing walls is "special" because of the peculiar critical 
point shift (||) which is somewhat anomalous and has at- 
tracted a lot of interest in recent years ||,^,|U|j2^] . In 
our opinion this point of view needs to be reconsidered; 
the mechanism of critical point shift becomes clear only 
for nonzero gravity where we find a scaling relation (j(^) 
completely analogous to that expected for the capillary 
critical point (|l|). 

Although they are currently restricted to two dimen- 
sions, DMRG techniques provide accurate results for 
studying equilibrium properties of large systems. The 
fact that the DMRG accuracy is the best for open bound- 
ary conditions [12 13 and that transfer matrices describe 



strips that are infinite along the transfer direction make 
the method best suited to study systems in contact with 



walls or with free surfaces. 
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